#PBAR_Zspec_DiscriminationPlots.py
#
#  Plots related to material discrimination
#
#  This code was moved over from PBAR_Zspec_ExamineData.py
#
#  5/4/2013, John Kwong

#import csv
import numpy as np
#import numpy.matlib
#import datetime
#import time
import PBAR_Zspec
#from pandas import Series, DataFrame
#from mpl_toolkits.mplot3d import Axes3D
#import pandas as pd
#from sklearn.lda import LDA
#from sklearn import svm
#import mlpy

#  SCATTER PLOT OF DISCRIMINANTS
fig = plt.figure()
detectorList = goodDetectorsList.copy()
##detectorList = detectorList[(detectorList>50) & (detectorList < 80)]
detectorList = detectorList[(detectorList>40) & (detectorList < 120)]

statsList = stats.keys()
# multibins
statName1 = 'transmission'
##statName2 = 'binSTD_binMean'
statName2 = 'binMean_g1'
##statName2 = 'binMean'
##statName2 = 'q_range_ratio'
##statName2 = 'multibin_10_ratio'
##statName2 = 'multibin_20_ratio_g'

disc1 = stats[statName1]
disc2 = stats[statName2]

# everything
cuttList = datasetGroupsIndices['Pb'], datasetGroupsIndices['Fe'], datasetGroupsIndices['Al'], \
           datasetGroupsIndices['Pb3Fe'], datasetGroupsIndices['Pb3Al'], \
           datasetGroupsIndices['Pb4Fe'], datasetGroupsIndices['Pb4Al'], datasetGroupsIndices['CC']

# everything except CC
cuttList = datasetGroupsIndices['Pb'], datasetGroupsIndices['Fe'], datasetGroupsIndices['Al'], \
           datasetGroupsIndices['Pb3Fe'], datasetGroupsIndices['Pb3Al'], \
           datasetGroupsIndices['Pb4Fe'], datasetGroupsIndices['Pb4Al']
cuttList = datasetGroupsIndices['Pb'], datasetGroupsIndices['Fe'], datasetGroupsIndices['Al'], \
           datasetGroupsIndices['CC']
cuttList = datasetGroupsIndices['Pb'], datasetGroupsIndices['Fe'], datasetGroupsIndices['Al']

# all lead
cuttList = datasetGroupsIndices['PbALL'], datasetGroupsIndices['Pb'], datasetGroupsIndices['Fe'], datasetGroupsIndices['Al']

cuttList = datasetGroupsIndices['PbALL'], datasetGroupsIndices['PbNOT']

##cuttList = datasetGroupsIndices['Pb'], datasetGroupsIndices['Fe'], datasetGroupsIndices['Al'], \
##           datasetGroupsIndices['PbLS'], datasetGroupsIndices['FeLS'], datasetGroupsIndices['AlLS']

for ii in range(len(cuttList)):  # cycle through list of materials
    cutt = cuttList[ii] # material cut

    temp1 = disc1[:,detectorList]
    temp2 = disc2[:,detectorList]

    kk = 0
    for jj in range(len(cutt)):  # cycle through list of datasets of that particular material type
        plt.scatter(temp1[cutt[kk],:], temp2[cutt[kk],:], marker = markerTypes[kk], color = plotColors[ii%7], s = (15 + 30*(ii/7)), label = datasetDescription[cutt[kk]])
        kk = kk + 1
        
plt.xlabel(statName1)
plt.ylabel(statName2)

# plt.yscale('log')
grid()
plt.legend(loc = 1, prop={'size':10})



#  SCATTER PLOT OF DISCRIMINANTS FOR THREE MATERIAls, cut in transmission
fig = plt.figure()
detectorList = goodDetectorsList.copy()
##detectorList = detectorList[(detectorList>50) & (detectorList < 80)]
detectorList = detectorList[(detectorList>40) & (detectorList < 120)]
statsList = stats.keys()
# multibins
statName1 = 'binMean_g'
##statName1 = 'multibin_0_10_g'
##statName1 = 'multibin_0_20_g'

##statName2 = 'binSTD_binMean'
##statName2 = 'binSTD_g'
##statName2 = 'q_range_g'
##statName2 = 'multibin_10_ratio_g'
statName2 = 'multibin_20_ratio_g'
##statName2 = 'multibin_0_20_g'
##statName2 = 'binMean_q50_g'

transmissionRange = np.array((10000, 15000))
transmissionRange = np.array((7500, 12500))
##transmissionRange = np.array((8500, 11000))
##transmissionRange = np.array((2200, 5200))

disc1 = stats[statName1]
disc2 = stats[statName2]
disc3 = stats['transmission']

cuttList = datasetGroupsIndices['Pb'], datasetGroupsIndices['Fe'], datasetGroupsIndices['Al'], \
           datasetGroupsIndices['Pb3Fe'], datasetGroupsIndices['Pb3Al'], \
           datasetGroupsIndices['Pb4Fe'], datasetGroupsIndices['Pb4Al']

cuttList = datasetGroupsIndices['Pb'], datasetGroupsIndices['Fe'], datasetGroupsIndices['Al']
# all lead
cuttList = datasetGroupsIndices['PbALL'], datasetGroupsIndices['Pb'], datasetGroupsIndices['Fe'], datasetGroupsIndices['Al']

cuttList = datasetGroupsIndices['PbALL'], datasetGroupsIndices['PbNOT']


for ii in range(len(cuttList)):  # cycle through list of materials
    cutt = cuttList[ii] # material cut
    temp = gainExtrapolated[0,88]/gainExtrapolated[0,:]

    temp1 = disc1[:,detectorList]
    temp2 = disc2[:,detectorList]
    temp3 = disc3[:,detectorList]
    
    kk = 0
    for jj in range(len(cutt)):  # cycle through list of datasets of that particular material type
        # plt.plot(temp1[cutt[kk],:], temp2[cutt[kk],:], marker = markerTypes[kk], markersize=12, label = datasetDescription[cutt[kk]], color = plotColors[ii])
        transCutt = (temp3[cutt[kk],:] > transmissionRange[0]) & (temp3[cutt[kk],:] < transmissionRange[1])
        plt.scatter(temp1[cutt[kk],transCutt], temp2[cutt[kk],transCutt], marker = markerTypes[kk], color = plotColors[ii%7], s = (15 + 30.0*(ii/7)), label = datasetDescription[cutt[kk]])
        plt.scatter(temp1[cutt[kk],transCutt], temp2[cutt[kk],transCutt], marker = markerTypes[kk], color = plotColors[ii%7], s = (15 + 30.0*(ii/7)), label = datasetDescription[cutt[kk]])
        kk = kk + 1

plt.title('Transmission range: (' + str(transmissionRange[0]) + ', ' + str(transmissionRange[1]) + ')')
plt.xlabel(statName1)
plt.ylabel(statName2)
# plt.yscale('log')
grid()


plt.legend(loc = 1, prop={'size':10})
# plot transmission 


#  SCATTER PLOT OF DISCRIMINANTS FOR THREE MATERIAls, cut in transmission, using collapsed array
fig = plt.figure()
correctx = 0
correctGain = 0
detectorList = goodDetectorsList[(goodDetectorsList>50) & (goodDetectorsList < 80)]

statsList = stats.keys()
# multibins
statName1 = 'binMean_g'
##statName1 = 'binSTD_g'
##statName1 = 'multibin_0_10_g'
##statName1 = 'multibin_0_20_g'
##statName1 = 'binMean_q50_g'

##statName2 = 'binSTD_binMean'
##statName2 = 'binSTD_g'
##statName2 = 'q_range_g'
##statName2 = 'multibin_10_ratio_g'
statName2 = 'multibin_20_ratio_g'
##statName2 = 'multibin_0_20_g'
##statName2 = 'binMean_q50_g'

transmissionRange = np.array((10000, 15000))
transmissionRange = np.array((7500, 12500))
transmissionRange = np.array((8500, 13000))
transmissionRange = np.array((2200, 5200))

datasetGroupList = np.array(('PbALL', 'Pb', 'Fe', 'Al'))
datasetGroupList = np.array(('PbALL', 'PbNOT'))

kk = 0
for ii in range(len(datasetGroupList)):  # cycle through list of materials
    group = datasetGroupList[ii] # material cut
    transCutt = (statsCollapsed[group]['transmission'] > transmissionRange[0]) & (statsCollapsed[group]['transmission'] < transmissionRange[1])
    plt.scatter(statsCollapsed[group][statName1][transCutt], statsCollapsed[group][statName2][transCutt], \
                marker = markerTypes[kk], color = plotColors[ii%7], s = (15 + 30.0*(ii/7)), label = group)
    kk = kk + 1

plt.title('Transmission range: (' + str(transmissionRange[0]) + ', ' + str(transmissionRange[1]) + ')')
plt.xlabel(statName1)
plt.ylabel(statName2)
# plt.yscale('log')
grid()

plt.legend(loc = 1, prop={'size':10})
plt.show()


#  3D SCATTER PLOT OF DISCRIMINANTS Using collapsed array
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
correctx = 0
correctGain = 0

detectorList = goodDetectorsList[(goodDetectorsList>50) & (goodDetectorsList < 100)]

statsList = stats.keys()
# multibins
statName1 = 'binMean_g'
##statName1 = 'binSTD_g'
##statName1 = 'multibin_0_10_g'
##statName1 = 'multibin_0_20_g'
##statName1 = 'binMean_q50_g'

##statName2 = 'binSTD_binMean'
##statName2 = 'binSTD_g'
##statName2 = 'q_range_g'
##statName2 = 'multibin_10_ratio_g'
#statName2 = 'multibin_20_ratio_g'
##statName2 = 'multibin_0_20_g'
statName2 = 'binMean_q50_g'

##statName3 = 'binSTD_g'
##statName3 = 'multibin_20_256_g'
#statName3 = 'multibin_0_20_g'
statName3 = 'multibin_20_ratio_g'

transmissionRange = np.array((10000, 15000))
transmissionRange = np.array((7500, 12500))
transmissionRange = np.array((8500, 13000))
##transmissionRange = np.array((2200, 5200))

datasetGroupList = np.array(('PbALL', 'Pb', 'Fe', 'Al'))
datasetGroupList = np.array(('PbALL', 'PbNOT'))

kk = 0
for ii in range(len(datasetGroupList)):  # cycle through list of materials
    group = datasetGroupList[ii] # material cut
    # plt.plot(temp1[cutt[kk],:], temp2[cutt[kk],:], marker = markerTypes[kk], markersize=12, label = datasetDescription[cutt[kk]], color = plotColors[ii])
    transCutt = (statsCollapsed[group]['transmission'] > transmissionRange[0]) & (statsCollapsed[group]['transmission'] < transmissionRange[1])
##    plt.scatter(statsCollapsed[group][statName1][transCutt], statsCollapsed[group][statName2][transCutt], \
##                marker = markerTypes[kk], color = plotColors[ii%7], s = (15 + 30.0*(ii/7)), label = group)
    ax.scatter(statsCollapsed[group][statName1][transCutt], statsCollapsed[group][statName2][transCutt], statsCollapsed[group][statName3][transCutt], \
                    marker = markerTypes[kk], color = plotColors[ii%7], s = (15 + 30.0*(ii/7)), label = group)
    kk= kk + 1

plt.title('Transmission range: (' + str(transmissionRange[0]) + ', ' + str(transmissionRange[1]) + ')')
plt.xlabel(statName1)
plt.ylabel(statName2)
ax.set_zlabel(statName3) 
# plt.yscale('log')

plt.legend(loc = 1, prop={'size':10})


# MATRIX SCATTER PLOT OF DISCRIMINANTS
statsCollapseList = np.array((
            'binMean_g1', 'binMean_q50_g1', \
            'multibin_0_20_g1', 'multibin_20_256_g1', \
            'multibin_20_ratio_g1', 'transmission'))
numberStats = len(statsCollapseList) - 1

cutt = statsMatrix[:,-1] == 0
statsDataFrame = DataFrame(statsMatrix[:, 0:5], columns = statsCollapseList[0:-1])


#ax = pd.tools.plotting.scatter_matrix(statsDataFrame.ix[cutt], alpha=0.2, figsize=(5, 5), diagonal='kde')
ax = pd.tools.plotting.scatter_matrix(statsDataFrame.ix[cutt], alpha=0.2, figsize=(numberStats, numberStats), diagonal = 'hist')

subplotx = numberStats
subploty = numberStats

#f, ax = plt.subplots(subploty, subplotx, sharex='col', sharey='row')

cutt = statsMatrix[:,-1] == 1
for ii in np.arange(numberStats):
    for jj in np.arange(numberStats):
        if ii == jj:
            ax[ii][jj].hist(statsMatrix[cutt,ii], bins = 20)
        else:
            ax[ii][jj].plot(statsMatrix[cutt,jj], statsMatrix[cutt,ii], '.r', alpha = 0.2)

###################
##   LDA PLOTS   ##
###################

for num in xrange(1, 1+len(statsCollapseList)):
    print("number of features %d" % num)
    # create short list of all sets that have this number of features
    cut = numberFeatures == num
    indices = where(cut)[0]  # references the super set
    successes = meanSuccess[cut]
    ranking = argsort(successes)[-1::-1]  # in order from highest to lowest; in subset space
    #show the top 5 sets
    for i in xrange(5):
        listlist = comboList[indices[ranking[i]]]
        print statsCollapseList[listlist], meanSuccess[indices[ranking[i]]], wMeanAll[indices[ranking[i]]]


# Scatter Plots of all the discriminants

for featureName in statsCollapseListALL:
    plt.figure()
    plt.grid()
    i = 0
    cutt = targetsAllTransmission == i
    plt.plot(statsCollapsed['PbNOT']['count'], statsCollapsed['PbNOT'][featureName], '.k', color = plotColors[i], alpha = 0.4, markersize = 10)
    i = 1
    cutt = targetsAllTransmission == i
    plt.plot(statsCollapsed['PbALL']['count'], statsCollapsed['PbALL'][featureName], '.k', color = plotColors[i], alpha = 0.4, markersize = 10)

    plt.legend(("Not Pb", "Pb"))
    plt.xlabel("Count")
    plt.ylabel(featureName)


# Scatter Plots of all the discriminants

categoryList = ['PbNOTNOT', 'PbALLALL', 'PbLS', 'FeLS', 'AlLS']
categoryList = ['PbNOTNOT', 'Pb3Fe', 'Pb3Al', 'Pb4Fe', 'Pb4Al']
categoryList = ['PbNOTNOT', 'Pb3Fe', 'Pb3Al', 'Pb4Fe', 'Pb4Al']
# categoryList = ['PbNOTNOT', 'Pb', 'Fe', 'Al']
categoryList = ['Pb', 'Fe', 'Al', 'Pb3Fe', 'Pb3Al', 'Pb4Fe', 'Pb4Al']

categoryList = ['PbALL', 'Fe', 'Al']
categoryList = ['Pb', 'Fe', 'Al']


for featureName in statsCollapseListALL:
    plt.figure()
    plt.grid()
    for (categoryIndex, category) in enumerate(categoryList):
        plt.plot(statsCollapsed[category]['count'], statsCollapsed[category][featureName], '.k',\
                 color = plotColors[categoryIndex], alpha = 0.4, markersize = 10, label = category)
    plt.legend()
    plt.xlabel("Count")
    plt.ylabel(featureName)
##    plt.title('Including Low Stats')
    
# Plot of the Discriminant Versus Count, each figure corresponds to training at a particularly count range.
comboMatch = [0, 1, 2, 3]
comboMatch = [1, 2, 3]
for (countIndex, countRange) in enumerate(countRangesList):
    for (comboIndex, combo) in enumerate(comboList):
        if combo == comboMatch:
            print(statsCollapseListALL[combo])
            wMean = wMeanAllAll[countIndex][comboIndex]
            dist0 = np.sum(np.matlib.repmat(wMean, statsMatrix.shape[0], 1) * \
                        statsMatrix[:,combo], axis = 1)
            plt.figure()
            plt.grid()
            for i in xrange(2):
                cutt = statsMatrix[:,-1] == i
                plt.plot(countMatrix[cutt], dist0[cutt], '.k', color = plotColors[i], alpha = 0.4, markersize = 10)
            plt.xlabel('Count')
            plt.ylabel('Distance')
            plt.title('Count Range: %s' % str(countRange))


# Plot of the Discriminant Versus Count, each figure corresponds to training at a particularly count range.
# Including Low Stats
comboMatch = [0, 1, 2, 3]
comboMatch = [1, 2, 3]
for (countIndex, countRange) in enumerate(countRangesList):
    for (comboIndex, combo) in enumerate(comboList):
        if combo == comboMatch:
            print(statsCollapseListALL[combo])
            wMean = wMeanAllAll[countIndex][comboIndex]
            dist0 = np.sum(np.matlib.repmat(wMean, statsMatrixIncludingLowStats.shape[0], 1) * \
                        statsMatrixIncludingLowStats[:,combo], axis = 1)
            plt.figure()
            plt.grid()
            for i in xrange(2):
                cutt = statsMatrixIncludingLowStats[:,-1] == i
                plt.plot(countMatrixIncludingLowStats[cutt], dist0[cutt], '.k', color = plotColors[i], alpha = 0.4, markersize = 10)
            plt.xlabel('Count')
            plt.ylabel('Distance')
            plt.title('Count Range: %s' % str(countRange))



# HISTOGRAMS WITH INCREMENTS OF VARIABLES
f, axAll = plt.subplots(1, featuresAllTransmission.shape[1])
for jj in range(featuresAllTransmission.shape[1]):
    ii = jj + 1

    disSubset = sum(np.matlib.repmat(wMean[topStatsIndices[0:ii]], features.shape[0], 1) * features[:,topStatsIndices[0:ii]], axis = 1)

    axAll[jj].hist(disSubset, bins = 30, label = 'distance subset')
    axAll[jj].set_title(str(jj))
    axAll[jj].set_title('Number of stats: ' + str(ii))
plt.legend()


# Maunally train on whole set

ldac = mlpy.LDAC()
cutt = countMatrix > 100
ldac.learn(statsMatrix[cutt,0:-1], statsMatrix[cutt,-1])
pred = ldac.pred(statsMatrix[:,0:-1])
sum(pred == statsMatrix[:,-1])
w = ldac.w()
dist0 = np.sum(np.matlib.repmat(w, statsMatrix.shape[0], 1) * \
    statsMatrix[:,0:-1], axis = 1)

plt.figure()
plt.grid()
cutt = pred.astype(bool)
plt.plot(countMatrix[cutt], dist0[cutt], '.k', alpha = 0.4, markersize = 10)

cutt = ~pred.astype(bool)
plt.plot(countMatrix[cutt], dist0[cutt], '.r', alpha = 0.4, markersize = 10)


plt.figure()
plt.grid()
cutt = statsMatrix[:,-1].astype(bool)
plt.plot(countMatrix[cutt], dist0[cutt], '.k', alpha = 0.4, markersize = 10, label =  'Not Pb')

cutt = ~statsMatrix[:,-1].astype(bool)
plt.plot(countMatrix[cutt], dist0[cutt], '.r', alpha = 0.4, markersize = 10, label = 'Pb')

cutt = ( pred.astype(bool) != statsMatrix[:,-1].astype(bool) ) & statsMatrix[:,-1].astype(bool)
plt.plot(countMatrix[cutt], dist0[cutt], 'bx', markersize = 7, label =  'Pb and Misclassified')

cutt = ( pred.astype(bool) != statsMatrix[:,-1].astype(bool) ) & ~statsMatrix[:,-1].astype(bool)
plt.plot(countMatrix[cutt], dist0[cutt], 'gx', markersize = 7, linewidth = 2, label =  'Not Pb and Misclassified')

xlabel('Count')
ylabel('Dist')



# PLOT THE SUCCESS AS A FUNCTION OF TRAIN SET #
plt.figure()
plt.grid()
plt.plot(successRateMean, label = 'Success on training set')
plt.plot(successRate0Mean, label = 'Success on all of data')
plt.xlabel('Training Set #')
plt.ylabel('Success')
plt.legend()

plt.figure()
plt.grid()
plt.plot(numberFeatures, succesRateMean, '.k', label = 'Success on training set', markersize = 15, alpha = 0.25)
plt.xlabel('# Features')
plt.ylabel('Success')
plt.legend()

